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Abstract 

A study of the quantum triangular billiard requires consideration of a boundary value 
problem for the Green's function of the Laplacian on a trianglar domain. Our main re- 
sult is a reformulation of this problem in terms of coupled non-singular integral equations. 
A non-singular formulation, via Fredholm's theory, guarantees uniqueness and provides a 
mathematically firm foundation for both numerical and analytic studies. We compare and 
contrast our reformulation, based on the exact solution for the wedge, with the standard 
singular integral equations using numerical discretisation techniques. We consider in detail 
the (integrable) equilateral triangle and the Pythagorean 3-4-5 triangle. Our non-singular 
formulation produces results which are well behaved mathematically. In contrast, while 
resolving the eigenvalues very well, the standard approach displays various behaviours 
demonstrating the need for some sort of "renormalisation" . The non-singular formulation 
provides a mathematically firm basis for the generation and analysis of periodic orbit ex- 
pansions. We discuss their convergence paying particular emphasis to the computational 
effort required in comparision with Einstein-Brillouin-Keller quantisation and the stan- 
dard discretisation, which is analogous to the method of Bogomolny. We also discuss the 
generalisation of our technique to smooth, chaotic billiards. 



1 Introduction 

The main focus of this paper is on contemporary issues in quantum chaos. In particular the 
convergence problems of periodic orbit expansions, the computational effort needed in their 
generation, and their status as the generalisation of Einstein-Brillouin-Keller quantisation to 
classically chaotic systems. To begin we consider the place our problem has in a broader 
context. 

1.1 The need for a semiclassical limit. 

Many problems in mathematical physics reduce to the solution of boundary value problems for 
partial differential equations or infinite matrices: wave propagation, diffusion, lattice statistical 
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mechanics, stochastic problems, thermodynamic formahsm in classical mechanics, an so on. 
Often, the "observables" can be expressed in terms of the eigenvalues and eigenvectors of some 
operator. The spectral problem of operators is of wide interest in mathematics as well. 

Unfortunately, these problems can only be solved in integrable cases, in which an n-degree 
of freedom problem reduces to n 1-degree of freedom problems. By "solve" we mean an effec- 
tive, efficient, mathematically well behaved means of calculation which affords an opportunity 
for intuitive understanding. The only way to proceed is to recognise that the problem has to be 
divided in two: the lower energy, long wavelength, portion of the spectrum and the high energy, 
short wavelength, portition of the spectrum. The former regime is approximately reduced to 
the diagonalisation of a finite matrix. In the latter regime one uses semiclassical approxima- 
tions. The stategy then is to develop the two techniques sufficiently and piece together their 
results to obtain a global understanding of the spectrum. 

For practical purposes the only global information one often needs is Weyl's density of 
states. Other times one simply requires the ground and first few excited states. It is then 
important to know at what point the spectrum crosses over into the "semiclassical" regime. 
It appears that the semiclassical regime extends quite far down. In the case of the disk, using 
asymptotic approximations for the Bessel function, one finds that a semiclassical approach 
gives acceptable approximations even for the ground state. For special examples such as the 
equilateral triangle and other problems solved by "Bethe Ansatz" the semiclassical approach 
gives the exact answer. For cases where global, detailed, spectral information is needed semi- 
classical methods are indispensible. 

1.2 Semiclassical techniques 

For historical reasons much of the development of semiclassical techniques has been driven by 
practical issues in optics and quantum mechanics. Such motivations have lead to mathematical 
developments such as the WKB method. Most of our motivation will come from quantum 
mechanics but optics will play a role. 

Semiclassical techniques have long had an important role in quantum mechanics. In the 
early days of the old quantum theory (and in fact in some contemporary problems) they 
provided the only quantisation technique. Failing to quantise more complicated systems, no- 
tably the helium atom, they were replaced by quantum mechanics. There, through the WKB 
method, they reappeared as approximation techniques in 1-dimensional or separable problems. 

Multidimensional, non-integrable quantum problems are much more difficult and are usu- 
ally treated with integral equations or diagonalisation in some basis. While providing the- 
oretical results to compare with experiment, such methods are computationally intense and 
afford little opportunity for intuitive understanding. Semiclassical techniques were first widely 
applied to such problems in optics Q. They have also been applied to quantum chemistry 
and more general problems of partial differential equations . 

1.3 Semiclassical quantisation. 

Although it was not widely realised, prior to the late 60s semiclassical quantisation was 
restricted to classically integrable systems using Einstein-Brillouin-Keller (EBK) quantisa- 
tion 0, ^, . It was then that Gutzwiller developed a method of semiclassical quantisation 
based on the calculation of the Green's function via a Laplace transform of the time-dependent 
propagator (heat kernel). For the latter a semiclassical approximation could be derived via 
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Feynman's path integral regardless of the integrablity of the corresponding classical system. 
Gutzwiller used the resulting periodic orbit expansion to determine the spectrum of the Kepler 
problem, and later, an anisotropic version of it. 

Around the same time, Balian and Bloch [0] obtained a similar result while extending 
a classical result of Weyl for the density of states of the Laplacian. They used a multiple 
expansion technique for the Green's function due originally to Carlemann. Weyl's formula 
is only relevant if one wishes to calculate certain sums over the spectrum. Balian and Bloch 
were interested in refining the Weyl result to incorporate the "shell structure" of the spectrum. 
They appear to have been motivated by certain spectral sums arising in nuclear physics. 

In recent years, mostly driven by resurgent interest in classical chaos, there has been re- 
newed interest in this work. The prime motivation has been to rehabilitate the old quantum 
theory by providing a quantization rule, analogous to EBK quantisation, for non-integrable 
systems. It is well known that periodic orbit expansions are absolutely divergent for chaotic 
systems due to an exponential proliferation of periodic orbits. Despite this, they have had some 
empirical success as rules for "quantising chaos". This motivatives the oft-stated hope that 
the series are conditionally convergent. The well-known properties of conditionally convergent 
series with regards rearrangment then begs the question: In what order should the sum be 
carried out? 

In an effort to overcome this difficulty, the original periodic orbit expansions have been 
refined in many ways. These refinements tend to be motivated from results in other areas 
of mathematics and physics: the theory of the Riemann zeta function, topology of hyper- 
bolic flows in classical chaos, the boundary integral method for partial differential equations, 
and scattering theory. These refinements include: the Riemann-Siegel lookalike formula and 
pseudo-orbit expansions |^] , cycle expansions [0] , surface of section techniques [0] and inside- 
outside duality [jl^] . 

Some refinement of Weyl's formulae in the calculation of spectral sums has been motivated 
from the Selberg trace formula for geodesic motion on hyperbolic manifolds [^, 13|. 

These techniques have been tested on a number of chaotic systems including: the hyperbola 
billiard ||l^, the helium atom the wedge billiard in a gravitational field 17 1, the Sinai 
billiard |18, the 3 disk system and smooth Hamiltonian systems |2l|]. Much 

progress has been made, as can most readily be seen in and where various modifications 
are compared. 

Despite the relative success of these improvements, they have all come from well motivated, 
but largely heuristic arguments and derivations. The well known mathematical difficulties are 
not directly tackled. An exception to this is some recent and very important work of Georgeot 
and Prange [^2[. Starting with the integral equations of Balian and Bloch they apply the 
Fredholm theory of integral equations to show, amongst other things, that the pseudo- 
orbit expansions for smooth billiards converge. The proper order of summation, of what is 
clearly a conditionally convergent series, is dictated by the Fredholm theory. Thus, at least for 
smooth billiards, it seems that the long standing questions regarding convergence have been 
settled, (well, maybe. See section \l.2\ ) 

Apart from convergence issues, there is the equally serious question of computational effort. 
It is fairly clear from the literature that, unlike EBK quantisation, periodic orbit expansions 
become more and more computationally intensive as one seeks higher and higher eigenval- 
ues. Furthermore this remains the case even with refinements which improve convergence p^ . 
While it is clear that exponential proliferation is the cause, it contradicts the correspondance 
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principle which suggests that the semiclassical approximation should improve with "larger 
quantum numbers". It may be, however, that this expectation is intertwined with our "inte- 
grable" intuition and needs to be reassessed in these "chaotic" times. 

Thus, despite the progress made, mystery continues to surround semiclassical quantisation 
of non-integrable systems. 



1.4 A role for the triangular billiard 

Our aim is to critically examine of the nature of semiclassical quantisation rules. To do this 
we seek out a model problem which is sufficiently simple to enable a first principles derivation 
which maintains reasonable standards of mathematical rigour. Whereas we do not seek, and 
are not capable of, a rigorous derivation, we aim at one which can plausibly form the basis for 
a rigorous approach. Our choice is the quantum triangular billiard. 

The quantum triangular billiard is well suited to a first-principles derivation. Being a 
billiard, the boundary value problem for the Green's function can be converted to an integral 
equation. This permits the application of a substantial body of standard theory. As with 
Georgeot and Prange, our main tool will be the Fredholm theory of integral equations. Most 
importantly, piecewise linear boundaries afford us much needed flexibility. Our discussion is 
immediately applicable to polygonal billiards. Although our notation will have an eye toward 
this generalisation, we will limit ourselves to triangles so as not be be distracted from our aim. 
The prospects for handling more general billiards, in particular those whose boundaries have 
smooth components, will be discussed in Section |7.2| . 

The quantum triangular billiard has been studied in a number of contexts. In the context 
of quantum chaology it has not attracted much attention because polygonal billards are not 
chaotic ||2^ (see however |26|]) so that one does not have exponential proliferation of periodic 
orbits. Despite this, it is still (generically) a non-integrable system, intermediate between 
integrable systems and chaotic ones, and so still of interest in its own right. It has also 



been considered in the contexts of eigenvalue degeneracies [p7|], energy level statistics |28], 
and modifications of periodic orbit expansions due to non-isolated and diffractive periodic 
orbits p^. Integrable triangles have eigenfunctions which are exactly expressed in terms of a 



finite superposition of plane waves via the Bethe Ansatz [30|. The general triangle has also 
been studied in an effort to extend the Bethe ansatz [^]. The Weyl formula for polygonal 
billiards has long been known |^2|. More recently this result has been extended to calculate 
determinant of the Laplacian 



1.5 Statement of problem 

To be more precise, the task we have set ourselves is the calculation of the spectrum of the 
Dirichlet Laplacian on a triangular domain. All spectral information conveniently resides in 
the Green's function which is the unique solution of the boundary value problem 

V'^T{x,Xi) + zT{x,Xi) = 0, xGA/{xj} (1.1a) 
lim </ VT{x,Xi) .dn = 1 (1.1b) 

T{x,Xi) = forxG^A (1.1c) 
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is a circular contour of radius e taken counterclockwise around the "source" point, Xj. A 
denotes the triangular domain and dA denotes its boundary. 

The unit source condition, (|l.lb| ), is usually specified as a Dirac 5-function inhomogeneity 
in ( |l.la ). We take the opportunity to remind the reader of this more precise description, even 



if it is somewhat old fashioned |34]. 



We are now in a position to be precise about the semiclassical limit. In this problem there 
are two natural length scales: the "de Broglie" wavelength \l \fz and any length L whose order 
of magnitude reflects the geometric distances over which the domain varies. For A it is the 
length of the sides. More complicated domains may provide other length scales which may or 
may not be a similar order of magnitude. The semiclassical limit is where the "extent" of the 
triangle is much larger than the "de Broglie" wavelength, \fzL ^ 1. Cast in this mathematical 
language, the semiclassical limit becomes a problem of asymptotic analysis. 

1.6 Standard integral equation approach 

For generic triangles the classical problem is non-integrable so that EBK is not applicable. 
The only way to obtain further analytic insight is to convert the boundary value problem 
into an integral equation using Greens' theorem. This is a standard technique in mathematical 



physics with a long history. It is widely used in electrostatics |34], diffraction and quantum 
scattering theory [^, both as a way to derive integral equations and also as a starting point 
for analytic and numerical approximations. It's (indirect) use in the study of quantum billiards 
began with Balian and Bloch We shall refer to this technique as the "standard approach" 
or the "free-based method" . 

For the Helmholtz equation various types of integral equations can be obtained. For the 
eigenfunctions one obtains an homogenous integral equation of the second kind. For the Greens' 
function one obtains an inhomogenous integral equation the second kind. Whilst there is a 
close connection between the two, the Fredholm theory gives the results needed to generate 
periodic orbit expansions only in the inhomogenous case. One is thus forced to study the 
spectral problem via the Green's Function. 

It is well known [^] (though sometimes overlooked) that in contrast with smooth billiards, 
the standard integral equations for polygonal billiards are singular. By this we mean that the 
requirements on the kernel demanded by Fredholm theory are not satisfied. Singular integral 
equations are known to be pathological. For example, the Lippmann-Schwinger equations 



for the quantum three body problem, which are singular, do not have unique solutions |3q] . 
In the literature on the Helmholtz equation, singular equations are known to yield "spurious 
solutions" [p7|, 37|. Although useful information can be extracted from these equations, the 



presence of pathologies creates the need for ad-hoc techniques, caution, and ultimately, doubt. 

The main result of this paper is the derivation of non-singular integral equations for the 
Green's function of the quantum triangular billiard. The derivation of non-singular integral 
equations from singular ones has a precedent in the work of Faddeev on the quantum three 



body problem [36|. The basic idea of Faddeev's work is that it is neccessary to treat part of 



the interaction in an explicit non-perturbative manner. Faddeev does this by eliminating the 



potential in favour of the T-matrix. As we will discuss in Section |2.8| a direct application of 
Faddeev's technique is not convenient for our problem. Instead we follow the spirit by basing 
a derivation on the Green's function for the wedge. 

Whereas most treatments of the quantum billiard problem depict the derivation of the 
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integral equations as straightforward and concentrate on the generation of periodic orbit ex- 
pansions, we shah show that for polygonal billiards this is not the case. Although we shall 
present a relatively compact and well motivated derivation, this does not reflect the tortuous 
path we took to reach it. 



1.7 Outline of paper 

We begin in Section ^ with an outline of the standard approach, presenting more detail than 
is customary in order to demonstrate some subtleties which are usually overlooked. The 
derivation will be accompanied by an elaboration of the relatively well known intepretation 
in terms of a "quantum Poincare section". We then demonstrate the singular nature of the 
standard integral equations and discuss ways in which they can be manipulated to obtain 
non-singular equations. We also emphasise the perturbative nature of the whole Green's 
theorem-based approach. We then warm up for the wedge-based approach by rederiving the 
standard integral equations from a perturbation about the half-plane Green's function. A 
similar derivation has recently been presented by Li and Robnik |39]. 



In Section y we outline the derivation of non-singular integral equations based on the exact 
solution for the wedge Green's function. The implementation of the non-singular formulation 
requires a detailed analytic and computational understanding of this solution. In Section ^ 
we give the necessary technical detail for the wedge-based kernel. In the appendix we give an 
outline of the derivation of the wedge Green's function. In Section ^ we present and compare, 
for purposes of validation, the results obtained from a standard numerical treatment of both 
the singular and non-singular integral equations. In particular we present the eigenvalues 
for the equilateral and Pythagorean 3-4-5 triangles by examining the zeros of the Predholm 
determinant. We show that the Predholm determinant for the discretised singular equation 
does not converge although it resolves the eigenvalues well. In contrast, the non-singular 
equation behaves well, as the Predholm theory demands. 

Given a suitable integral equation, the derivation of periodic orbit expansions is more or 
less mapped out by Georgeot and Prange. Let us briefly outline their discussion. The Predholm 
theory expresses the solution of an integral equation in terms of a resolvent kernel Hx and a 
Predholm determinant A(A) in the form 



1 



b 



xit) = xoit) + I Hx{t, s)xois) ds (1.2) 



|23|. The eigenvalues are then the zeros of the Predholm determinant, which one identifies 
as a zeta function. The Predholm determinant can be expanded in series using the identity 
logdet(/ — K) = Trlog(l — K) [^8|. The expansion of log(I — K) in powers of K has a finite 
radius of convergence. However exponentiating the series and re-expanding gives the standard 



Predholm series which has an infinite radius of convergence [23, Theorem 6.5.2]. A semiclassical 
expansion of the traces Tr:{K^), which appear in this series, then provides a convergent periodic 
orbit expansion for this zeta function which one identifies with the psuedo-orbit expansion. 

At this point it is necessary to make a terminological aside. The term "periodic orbit 
expansion" can refer to the original Gutzwiller series or any of its refinements. The context 
of discussion allows one to use it without qualification. However, in order to emphasise the 
firm foundation provided by the Predholm theory, we feel that it is perhaps appropriate to 
occasionally use the term "semiclassical Predholm series" . 
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In Section |^ we use the Predholm theory to analyse the convergence of resulting periodic 
orbit expansion. Here it will become clear that while Fredholm theory solves the convergence 
problem the number of terms needed for convergence grows like \fE. The semiclassical ap- 
proximation of the n-th term in the series expansion of the Fredholm determinant requires 
all ra-bounce periodic orbits. For the triangle the number of such orbits grows polynomially 
in n. The computation effort in the periodic orbit expansion is then polynomial in \fE. In 
contrast the effort for chaotic systems is exponential. We compare this effort with the O(n^) 
effort needed in the numerical discretisation of the integral equations. This comparision clearly 
defines how far we are from the goal of emulating EBK quantisation for chaotic systems. 

We conclude with a discussion in Section |^. Here we present the singular equation as a 
simple model system which requires "renormalisation" . We note the presence of logarithmic 
singularities in the Fredholm series making it comparable to perturbation series in quantum 
field theory. The wedge based theory is then interpreted as a first principles, mathematically 
firm, "renormalisation scheme" where any ac?-/ioc "regularisation" and "subtraction of infini- 
ties" is absent. We then present arguments that cast doubt on the validity of the Balian-Bloch 
integral equations for smooth billiards. We then suggest ways in which our construction can 
be carried over in principle to smooth billiards. We also discuss the relationship between the 
Fredholm series approach and other modifications of Gutzwiller's original trace formula. We 
also present our construction as a way to derive, in a mathematically satisfactory way, Keller's 
geometrical theory of diffraction Q. 



2 Standard Approach: Singular Integral Equations 

2.1 Standard pre— integral equation 

The standard approach is based on the free Greens function. This is the solution of the 
Helmholtz equation on the Euclidean plane which decreases at infinity and has a unit source. 

Go{x,Xf) = -'-H!,'\^\x-Xf\) (2.1) 



where H^\z) is the zeroth order Hankel function of the first kind [40|. 

Our main point is that this is a choice and that for polygonal billiards it is a poor choice. 
Let us thus begin with any Greens' function, i.e. any function F{x,Xf) which is a solution of 
the Helmholtz equation with a unit source. The domain on which F{x,Xf) is defined and the 
boundary conditions it satisfies are totally arbitrary. 

Applying Green's second identity, 

((/)V V - V'VV) dx = / {4>V4^ -i{jV4>) .dn, (2.2) 

D JdD 

with (p = T{x,Xi) and ip = F{x,Xf) to the domain D given in Figure ||, taking the width of 
the channels 5 to zero, and using continuity in the usual way, one obtains 

(/ -/ -/ ][T{x,Xi)VF{x,Xf)- F{x,Xf)VT{x,Xi)].dn = 0. (2.3) 
\JdA Jc^^ Jc^^. J 

and are circular contours of radius e. 
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Figure 1: Contour for derivation of pre~integral equation. For convenience distances t are 
defined with respect the nearest vertex in a counterclockwise sense. The heads of vectors are 
denoted using points. The integration region D is shaded and the sense of dD is shown. 



In the integral over Cx- all functions, except VT, are continuous and can be approximated 
to 0(e) by their values at Xi. Using this, the unit source condition, and taking e — > we obtain 



lim (t [T{x,Xi)VF{x,Xf) - F{x,Xf)VT{x,Xi)] .dn 



-F{xi,Xf) 
-F{xi,Xf) 



lim / \7T{x,Xi) . dn 



(2.4) 



The integral over Cxf is evaluated similarly. 

Using the Dirchlet boundary condition satisfied by T{x, Xi) we finally obtain 



T{xf,Xi)=¥{xi,Xf)- / F{x,Xf)VT{x,Xi) .dn 

JdA 



(2.5) 



The integration in ( p.SP is concentrated on the boundary of the domain but also involves values 
of the unknown function at arbitrary points in the interior. Interpreted as an integral equation 
on A it has a very singular kernel. The standard approach instead regards this as an equation 
which determines the solution in the interior given it's normal gradient on the boundary. For 
this reason we refer to ( |2.5| ) as the pre-integral equation. Let us pause for a moment to interpret 
this equation. 



2.2 Poincare section interpretation 



The boundary of the domain in (2^) plays a distiguished role very similar to its role in the 
classical problem. In the classical problem Birkhoff defined a Poincare section for billiards 
consisting of pairs {q,p) giving the position and tangential momentum at each collision with 
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the boundary. A particular classical orbit is then represented in the Poincare section by a 
sequence of points, {qmPn),^ = 0,1,2... ,00. Knowing that the motion between collisions 
is free we can explicitly give the full orbit from this sequence. For billiards the first return 
map, Qn+i = fiQn,Pn),Pn+i = giQn,Pn), be obtained explicitly. This "surface of section" 
construction then reduces the problem to a simpler dynamics, that of a map, on a space of 
reduced dimensionality. 

This construction is identical in spirit to the above interpretation of the pre-integral equa- 
tion. Rather than consider the space of all functions defined on A as dynamical variables we 
have reduced the dimensionality of the unknowns to the space of functions on dA. All we need 
an integral equation on this space to play the role of a "quantum first return map" and whose 
unique solution is the actual boundary normal derivative. 

This analogy is of course not direct because of the fundamental differences between classical 
and quantum mechanics. In quantum mechanics there is no sensible concept of "phase space". 
The state of the system, in the Schroedinger form of quantum mechanics, involves only the 
interior and boundary configuration coordinates. Furthermore q is the independent variable of 
the wave function rather than a dynamical variable q{t). 

It would appear that a "quantum first return map" can be obtained in a straightforward 
manner by taking the gradient of ( |2.5D with respect to Xf^ evaluating the limit as Xf tends to 
the boundary and projecting out the normal derivative. We shall refer to this process as the 
boundary limit. 

2.3 The boundary limit 

The boundary limit is in fact rather subtle. Let us now take up the standard approach in 
which the free Green's function is chosen for F(x,Xf). 

At this point it is necessary to introduce some notation. Refering to Figure ^ we locate 
the verticies of the triangle at Pq,Pi and P2- The sides of the triangle are labelled with the 
integers 0,1,2 so that 9 A = + 1 + 2. The length of side is 



For convenience we parameterise each segment of the boundary using the distance from 



^0 1 




Figure 2: Notation: Free-based integral equation 
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the vertex immediately clockwise: 



\P1 


-Po\ 




-Pi 


\P2 


-Pl\ 


Po 


-P2 


\Po 


-P2\ 



Then we define 



ro{s) = Po + 
ri(s) = Pi + 
r-2(s) = P2 + 



X^{t) = VT{r^{t),Xi).n^ 
Xf,{t) = 2VGo{r^{t),x,).n^ 
K^^{t,s) = -2VGo{r^{s),r^{t)).n^ 



0< s<lo 
< s < /i 
< s < /2 

< t < / 
< t < / 



<t < l^i, < s < 



(2.6a) 
(2.6b) 
(2.6c) 



(2.7a) 
(2.7b) 
(2.7c) 



where /U = 0, 1, 2, = 0, 1, 2. The need for the mysterious factors of 2 will soon be divined. 
Formally taking the boundary limit towards side 0, we obtain 



2X0 (t) = X^{t) + Kooit, s)Xois) ds 
Jo 



+ / KQi{t,s)Xi{s)ds+ / Ko2{t,s)X2{s)ds 







(2. 







In order to examine the boundary limit in detail we first present a diagramatic notation for 
(|2.8| ) using the 'Teynman rules" in Figure |^. Diagrammatic notation is totally equivalent to 
the algebraic notation and provides a convenient tool for generating perturbation expansions. 
In the present context it clearly elaborates the structure of the integral equations. It will also 
prove very useful and natural in the generation of periodic orbit expansions. 






Figure 3: Feynman Rules for translating between integral equations and diagrammatic rep- 
resentations. The arrows give the "flow" from Xj to in the integral equations. They are 
needed to distinguish between Kfj_^{t,s) and Ky^{t,s). 



The translation between notations is given by the usual rules 
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• Associate each factor in the integrand with an edge as given in Figure ^. 

• The configuration space dimensions in the arguments of the integrand functions can be 
labelled on the triangle as shown. However these are fairly obvious and usually omitted. 

• An integration along a side is associated with each vertex where two edges meet. 

In Figure ^ we give a graphical representation of ( |2.8| ) . The application of the Feynman 
rules is shown for part of the equation. 



2 





Figure 4: Graphical representation integral equation in naive boundary limit. Also shown is 
an example of the application of the Feynman rules. 

At this stage we require the explicit form of the integral kernel 

K,At,s) = -'^Hi'\V-z\ A,. (M) I ) "^f^^'^'l ^^^^ (2.9) 

Here A.^y{t, s) = r^{t) — r^{s) and h[^^ (z) is a first order Hankel function of the first kind |4C]. 
The subtlety in the boundary limit arises because 

Apart from this the kernel is well behaved. The singularity occurs when |A^,y(t, s)| is 
much smaller than the de Broglie wavelength so that it appears in the quantum regime. This 
distance is explicitly represented in the diagrammatic notation by the length of the segment 
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corresponding to the "free propagator", It is thus easy to see when a singularity is 

encountered in the integration. 

Let us now go back and do the boundary Umit towards side in a more careful manner. In 
boundary limit we take aij — > ro{t). Let us make the limit explicit by letting xj = 'ro(t) — eno. 
We can readily "deform" our diagrammatic notation to accomodate this. In particular, in 
Figure Q we can denote this deformation by simply displacing the final point slightly above the 
boundary. 

It is easy to see that the boundary limits of the left hand side and the inhomogenous 
term exist by continuity, since Xi is strictly in the interior. It is also clear that, provided 
we stay away from the corners to keep the distance of free propagation bounded from below, 
the integrations over boundaries 1 and 2 have uniformly bounded integrands. They are then 
uniformly convergent and the boundary limit can be taken inside the integrals to give the 
formal result. The same is not true for the integration over boundary 0. Here it is clear that 
the integrand develops a singularity in the boundary limit due to the presence of arbitrarily 
short distance propagation. The consequent lack of uniform convergence means that the formal 



boundary limit is not correct. This observation has already been made by Li and Robnik |3£]. 
Their remedy is to perturb about the half-plane. Their treatment fits into our general scheme 
and we shall describe it in Section |2.7| . 

The boundary limit towards side can still be done, but care is required. In explicit form, 
we need to calculate. 

rlo . 

-Xo{s)ds 



lim 



s)2 + e2) 



s)2 + e2 



(2.11) 



The integrand tends to zero as e — > for t ^ s. Because the singularity in ( |2.10| ) reflects the unit 
source condition satisfied by the free Green's function we expect to see (5-function behaviour. 
This can be revealed explicitly by scaling the dummy variable according to s = t + ex to obtain 



lim 



(lo-t)/e 
■t/e 



1 



Xo(t + ex) e dx 



(2.12) 



Expanding Xo(t + ex) as a Taylor series in e and using the small argument form for the Bessel 
function we can explicitly evaluate the integral for small e to obtain 



1 filo~t)/e 

lim — / 

27r y„t/. 



x^ + 1 



+ O(eloge) } dx 



Xo{t) = -Xo{t) 



(2.13) 



We see then that in the boundary limit the integrand behaves as ^6{x). 

Repeating the calculation for the boundary limit toward sides 1 and 2 (or using a cyclic 
permutation of labels) , we obtain 



Xoit) = X^{t) + r Koiit, s)Xiis) ds + Ko2it, s)X2is) ds 
Jo Jo 

rlo rh 



X2{t)=Xl{t)+ [" K2oit,s)Xo{s)ds+ r K2i{t,s)Xi{s)ds 
Jo Jo 



(2.14a) 



Xi{t) = X^{t) + [ " Kioit, s)Xo{s) ds+ r Kuit, s)X2{s) ds (2.14b) 
Jo Jo 

(2.14c) 
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This is a coupled form of the standard integral equation. 

From (|2.9|) we observe that KQo{t,s) = because ro{t) — ro(s) _L no- It would be quite 
easy, in a first treatment, to ignore the fact that this is strictly only true for t ^ s and use it in 
the formal boundary limit. This would lead to the same set of integral equations, except that 
the factors of 2 in (2.7) would be absent. We made this mistake in our early investigations and 
are quite confident that in doing so we have tread a well worn path. 



2.4 Other treatments of the boundary hmit 

Balian and Bloch 0[I] derived the standard integral equations using potential theory. Their 
pre-integral equation (11.10) with (11.12) is of a different form to ours due to the choice of 
a double layer potential. We cannot understand how such an object can be obtained from 
Green's theorem. The multiple scattering series which they give is nevertheless identical to 
the one obtained by iterating our equations. Potential theory is essentially the application 
of Green's theorem to electrostatics. Because the Poisson equation is very different to the 
Helmholtz equation we feel that, while potential theory may give correct results, the delicacy 
of the problem demands the direct use of Green's theorem. 

Other authors p3, 27] derive the standard integral equations from the Helmholtz form of 



Greens' second identity. 



/ ds' {ij{r')n' .Vr'Go{r,r') - Goir,r')n' .Vr'i^ir')} 

JdD 



>(r) r G D/dD 



(2.15) 







r ^ D 



The mysterious factor of ^ can be understood as being due to the indentation of the contour 
with a semicircle when r appears on the boundary. 

The derivation proceeds by differentating ( p. 151 ). As our discussion shows, there is clearly 
some sort of boundary layer behaviour so that ( p.l5[) needs to be supplemented with a statement 
of this behaviour. We suspect that performing the differentiation within the boundary layer 



is justified. As such we feel that Li and Robnik |39| are not strictly correct in describing the 
approach as in error. As the reader now has our discussion of the boundary limit, as well as 
the approach or Li and Robnik, we feel no particular need to elaborate any further. 

Our study the boundary limit was motivated by the hope, forlorn as it turned out, that by 
doing the problem "carefully" the singular nature of the integral equations could be removed. 



2.5 Singular nature: Predholm theory 

In taking the boundary limit it is convenient to use a notation which identifies each side of the 
triangle explicitly. This coupled form is convenient practical calculations because it provides 
a symbolic dynamics for quantum progagation. In the semiclassical analysis this will reduce to 
the natural symbolic dynamics for the classical propagation. 

For theoretical purposes, in particular to apply the Fredholm theory, it is convenient to 
convert these coupled integral equations into a single integral equation. To do this we first 
define the cumulative arc length to each vertex from Pq: ctq = 0, cr^+i = + = 1, 2, 3. 
We can then define x{a) = r^[a — a^) and n{a) = for < a < cx^+i giving the position 
and normal vectors in the Birkhoff co-ordinate a. 



13 



Then we define 



X{a) = VT{xia),x^) . n{a) (2.16a) 
Xoia) = 2VGo{x{a),Xi).n{a) (2.16b) 
K{t, a) = -2VGo{x{T),x{a)) . n(T) (2.16c) 



and rewrite ( 2.14 ) in the more compact form 

X,it) = XO(t) + / K,At, s)X,is) ds (2.17) 



where we have redefined K^^{t, s) = 0. This plays the role of a quantum pruning rule, resulting 
from the boundary limit process. We then deduce 

K^^{t,s) = K{a^ + t,a^ + s) (2.18a) 
X^{t)=X{a^ + t) (2.18b) 
Xl{t)=X\a^ + t) (2.18c) 

( |2.17| ) can then be translated into an uncoupled form 

X{t) = X°(r) + / ' K{t, (j)X{a) da (2.19) 

J (TO 

In this form we can readily apply the Fredholm Theory. The important quantity here is 
the norm of the kernel. 



\\K\\ = Y^Tr(Ki^t) = J f ' da f ' dT\K{T,a)\^ (2.20) 
y J (TO J (70 

where ^ denotes the Hermitian adjoint. In order to estimate the norm we revert to the 



coupled form using ( p. 18]) 

Y^Y. ds\K^u{t,s)\^ (2.21) 

f,=Ou=0-^0 Jo 

The norm is the sum of 6 positive numbers. It is clear from the diagrammatic notation 
that all 6 terms contain short distance propagation. Let us now show that these short distance 
contributions cause the norm to diverge. Consider the case 

[" ds r dt\K2o{t,s)\'^ (2.22) 
Jo Jo 

The other 5 cases will follow similarly. Let t ^ I2 — t so that both distances are refered to 
vertex 0. 

The integration region consists of the rectangle (0, Zq) ® (0,/2)- In this rectangle draw a 
quarter circle of radius 6 <^ l/\/^ centered on the origin. This divides the integration into 
two domains. In the bounded domain outside the quarter circle the integrand is bounded and 
continuous so that the integral exists. Further subdivide the remaining quarter circle using 
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a smaller quarter circle of radius e centered on the origin. Consider the integration over the 
quarter annulus so defined. Because A <^ this region we can use ( 2.10| ). Using polar 

coordinates in the (t, s) plane one readily shows that the integration over the annular region is 



^/(0o)log- (2.23) 

where I{(po) is some unimportant integral. 

The norm diverges logarithmically and the Fredholm theory does not apply. In fact a portent 
of the singular nature appears in a more careful examination of the boundary limit. In ( p. 12 ) 
it was implicitly assumed that Iq — t and t were not 0(e). If we had taken the boundary limit 
whilst keeping t/e constant we would have obtained a different, e-dependent result. (This 
gives a quantitative measure of what we meant by the phrase "provided we stay away from 
the corners".) This involves the wave function a distance e into the interior of the domain 
so that resulting integral equation is defined on A rather than dA and is not the "quantum 
first return map" we seek. The reader may be able to convince herself that there is no direct 
remedy for this boundary layer behaviour. 

2.6 Perturbative character of boundary integral method 

Confronted with no apparent way to overcome this singular behaviour in a mathematically 
satisfactory manner, one is forced to reassess the derivation. The first thing the reader must 
appreciate is that the method presented is a perturbative method. 

In quantum and classical mechanics complicated Hamiltonians are treated by spliting off 
an "interaction" term V so that H = Hq + V and the problem corresponding to Hq can 
be explicitly solved. One then constructs a perturbation series based on this solution. The 
conversion of the boundary value problem to an integral equation presented above is directly 
analogous to this common technique. Here the simplification is not achieved by spliting the 
Helmholtz operator but by simplifying the boundary conditions. In an operator theoretic 
sense there is no distinction between the two simplifications. The result, however, is not a 
perturbation series but an integral equation so that it is a non-perturbative result from a 
method with perturbative character. 

The standard derivation obtains an "Hq" with a very brutal reduction of the boundary 
value problem. Our investigations have led us to conjecture that: 

It is in principle not possible to obtain the full triangle Green's function perturbing 
about the free Green's function 

If this conjecture is correct then the next step is obvious: 

One must perturb about a Green's function which reflects more closely the basic 
structure of the full solution 

At this point we require a catalog of exactly known Green's functions. 

• For domains which tile Euclidean space under reflection, the Greens function can be 
constructed by the method of images. 

• For the wedge Sommerfeld generalised the method of images to give the Greens function 
as a contour integral. 
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• For the circle and other separable problems one can obtain eigenfunction series for the 
Greens function. Semiclassical studies then require use of complex angular momentum 
techniques (aka. Sommerfeld- Watson transform, Regge poles pOt|). 

The reader familiar with these solutions is in a position to appreciate that the availability 
and relative simplicity of the method of images solutions and the wedge solution was the prime 
motivation in our decision to study the quantum triangular billiard. 



2.7 Perturbing about the half-plane 

The next simplest Green's function about which to perturb is the Dirichlet Green's function 
on the half-plane Dh 

H{x,Xf) = Goix,Xf) - Goix,xf), (2.24) 

obtained using the method of images, x^- is the mirror image of Xf in the boundary dDu- 

Li and Robnik [^] have presented a half-plane based derivation for smooth billiards and 
shown that the standard integral equations are recovered. Let us briefly present this derivation 
for the triangle both as a warmup to the wedge case as well as an elaboration of Li and Robnik's 
result. 

First extend side to infinity in both directions. Denote the half-plane Green's function 
which vanishes on this line by HQ(x,Xf) Using this Green's function in ( ^.5| ) gives 

T{xf,Xi) = Ho{xi,Xf)- [ Ho{x,Xf)VT{x,Xi).dn (2.25) 

Jl+2 

The Dirichlet boundary condition satisfied by HQ{x,Xf) eliminates the delicate integration 
over side 0. Because of this the boundary limit becomes straightforward, provided, once again, 
that we are not too close to the corners. Perturbing about the half-plane Green's functions 
which vanish on sides 1 and 2 and taking the boundary limit we recover the standard integral 



equations (2.14). The "mysterious" factors of 2 in (2.7) arise here because the image term 
gives a contribution to the kernel which is identical to that of the free term. 

Note that Sieber uses a perturbation about the quadrant Green's function in an application 



of the boundary integral method to the hyperbola billiard [|14|. 



2.8 Direct removal of singular nature 

It is possible to derive a non-singular integral equation from the singular integral equation. 
Riddell does this by separating out the singular part of the kernel and treating it exactly 
via an auxiliary integral equation. Riddells' technique is somewhat cumbersome and is not 
particularly intuitive. Faddeev ||3^, in his treatment of the quantum three body problem, used 
a more systematic version of this approach. 

We began our studies trying to following Faddeev's method. We were able to successfully 
adapt Faddeev's techniques to our present situation by splitting the (matrix) kernel into three 
pieces where each piece involved only collisons between two sides. This is similar to the 
decomposition of the potential energy in the three body problem as V12 + V23 + V31. One can 

^As this is the only place where the half-plane solution is used, no confusion with the notation for the Hankel 
functions should arise 
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then follow Faddeev fairly closely, defining T-matricies for each "wedge" of the decomposition 
and using these to derive non-singular integral equations. We found the result unsatisfactory 
because we were unable to obtain a convenient form for the T-matricies. 



3 Wedge— Based Approach: Non-Singular Integral Equations 

The next most natural way to proceed is to base a derivation on the Dirichlet Green's function 
of the wedge. 



3.1 Historical overview 

The scattering problem for the wedge is of considerable importance in optics where it is virtually 
the only non-trivial example in which Maxwell's equations can be solved exactly [^]. It was 
first considered by Sommerfeld who approached it via a generalisation of the method of images 
to Riemann surfaces. Sommerfeld's techniques were simplified and applied to other problems 
by Carslaw It is a classical problem in mathematical physics and is treated in many text 
books as well as in the research literature |42]. Attention is often given to the special case of 



the edge where the wedge angle is 27r. 

The wedge solution has been used as the basis for analysis in a number of studies. Kac 
used the solution of the diffusion equation on a wedge to determine the constant term in Weyl's 



formula for the density of states. More recently Aurell and Salmonson |33| have presented a 
similar argument in their study of the determinant of the Laplacian on triangular domains. 
The wedge solution has also appeared in the recent work of Pavloff and Schmit which 



incorporates wedge diffraction in periodic orbit expansions for the triangle. Diffraction effects 



in periodic orbit expansions were first considered by Vattay, Wirzba and Rosenqvist |43] using 
Keller's geometrical theory of diffraction . In these studies scattering off verticies involves a 
factor which originates in the asymptotic form of the scattering solution for the wedge . 

In our case we require the Greens function for the wedge satisfying a unit source condition 
rather the condition of an "outgoing wave" at infinity. This can be obtained using Sommerfeld's 
method [^] although it was originally obtained by Macdonald in 1902 using other methods |45[| . 
The solution is of same form as that of the scattering problem except that the plane wave is 
replaced by a Bessel function of a somewhat more complicated argument. Because the Green's 
function has been considered to a lesser degree we shall give some details. Despite appearances 
it is just as workable as the scattering solution. In the appendix shall give an overview of its 
derivation. 



3.2 Naive, wedge— based integral equations. 

Refering to Figure || let us label the wedge solutions according the vertex they share with the 
triangle. Perturbing about Wo{x,Xf) and using ( ^^ ) gives the pre-integral equation 

T{x,Xi) = Woix„Xf) - j^Wo{x,Xf)VT{x,x,) .dn (3.1) 

where the Dirichlet boundary condition satisfied by the wedge removes two integrations. 

One can now take the boundary limit to either side or side 2 without encountering the 
short distance singularity in the propagator. However as Figure ^ illustrates, taking the limit 
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too close to corners 1 or 2 produces short distance propagation. Because the wedge solution 
satisfies a unit source condition its short distance behaviour is the same as the free Greens 
function so that, once again, short distance propagation causes the norm of the resulting kernel 
to diverge. 



Figure 5: Corner divergences in the naive wedge-based derivation: perturbation about the 
wedge solution based on vertex 0. 

Carrying through the derivation regardless, one obtains 6 coupled equations for 3 unknown 
functions. This is an overspecified system, as is most readily seen by discretising the result 
to obtain a problem of linear algebra. It is possible to remove the short distance problems by 
taking advantage of this overspecification. 

3.3 "Regularising" the singular nature 

We begin by adding an extra vertex on each side of the triangle so that Figure |2| becomes 
Figure ^. In this manner we regard the triangle as a degenerate 6-sided polygon. The manner 
in which the sides are subdivided is arbitrary but, as we will discuss later, there will be an 
optimal subdivision. 

The boundary is now parameterised in 6 segments 

^/^(^) = n+ |p^^'"pV 0<'<^^ (3-2) 

where fi = 0,1,2,3,4,5 and Pq = Pq in a natural cyclic manner. We denote the outward 
normals by tiq = rii, n2 = and n4 = n4. 

The three wedge solutions about which one can perturb are now associated with verticies 
0, 2 and 4. Let us consider the perturbation about Wo{x,Xf). 

Define 

X^{t) = VT{r^{t),Xi) .n^ Q<t<l^ (3.3a) 

Xlit) = 2VW^oK(t), Xi).n^ 0<t<l^ (3.3b) 
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Po: oPi* 1 \P2 

lo h 
Figure 6: Notation: Wedge-based integral equation 

where fi = 0,1, 2, 3, 4, 5 and 

Kf,u{t, s) = -2VWo{r^is), r^{t)) .n^ 0<t<l^, 0< s < (3.4) 

where /U = 0, 1,4, 5 and u = 2,3. (|2.5| ) then gives the pre-integral equation. 

T{xf,x,) = Woixi,Xf)+ Wo{r2{s),Xf)X2is)ds+ [ ' Wo{r3{s),Xf) X^is) ds (3.5) 

Jo Jo 

There are now 4 possible boundary limits. Defining a diagramatic notation analogous to 
that given in Figure ^ we illustrate 2 typical boundary limits in Figure |^ Taking the limit to 
sides and 5 it is clear that there is a lower bound on the propagation distance. This is not the 
case taking the limit to sides 1 and 4. We take advantage of the overspecification, only using 
the non-singular limits. 

Perturbing about the wedge solutions associated with vertices 2 and 4 in an analogous 
manner gives a set of six coupled integral equations 

X^{t) = X0(t) + V / K^,{t, s)X^{s) ds (3.6) 
u=o-'^ 

Of the 36 possible K^y only the following 12 are non zero 

i^^,^+2,-?^At,M+3 M even (3.7a) 
K^^^+3,K^^^+4, /u odd (3.7b) 

where addition is modulo 6. These quantum pruning rules are easier to remember using the 
diagrammatic form of the integral equations. The kernel iC^jy is constructed using the wedge 
solution whose vertex is closest to the side corresponding to the first index. 

Perturbing about the wedge thus gives us a closed set of non-singular integral equations 
for the unknown normal derivative. The Predholm theory guarantees that these equations 
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^0 Pi 1 P2 Po ° Pi ' P2 

(a) Non-Singular boundary limit (b) Singular boundary limit 

Figure 7: Boundary limit of wedge pre~integral equation ( |3.5| ) to two typical segments. The 
integration over sides 2 and 3 is schematically given to show the absence (presence) of short 
distance propagation in the boundary limit to (1). 



define the unknown normal derivative uniquely. Substituting this solution into the pre-integral 
equation then gives the solution of the boundary value problem (1.1). 

In a sense, the need to start from the wedge solution is quite appealing because it gives 
a "bootstrap" character to the solution of boundary value problems: the triangle is solved in 
terms of the wedge problem; the wedge in terms of the free Green's function. One is led to 
think how the triangle solution might be used in constructing even more complicated solutions 
(perhaps the tetrahedron). 



3.4 Validation 

Although the derivation has been straightforward and the result is more or less clear, it is 
nevertheless important to validate our approach. The first step towards this is to carry out a 
standard numerical discretisation of the integral equations, calculating the Fredholm determi- 
nant and thus the spectrum. It is important to reproduce exactly known spectra and compare 
the results for general triangles with the singular approach. We shall do this in section |5[ 

In our numerical work we would like to be assured that any discrepencies, any problems, 
are entirely due to the formulation of the integral equations. It is then important to be able to 
calculate the kernel of the wedge-based approach with well-controlled approximations. This 
requires a detailed understanding of the wedge kernel. It is to this that we now turn our 
attention. 
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4 Wedge Kernel 



In appendix ^ we outline a derivation of the wedge Green's function and obtain the following 
form of the wedge kernel ( [A. 15 ) 



-VGw{rf,0;u,di)- (-§ 



1 



(4.1) 



where 



(prj jc 
Piei,cj);u;) 



1 



1 



1, 1 Z^'^/ , 
cosh a = - — H 

2 \ri Vf 



\{uj) = V cosh 



a — cos oj 



(4.2a) 

(4.2b) 

(4.2c) 
(4.2d) 



where the natural polar co-ordinate system has been used. h[^^ is the first order Hankel 



function of the 1st kind. The square root in (4.2c) is the principal square root with a branch 
cut on the negative real axis. The contour C together with the analytic structure of the 
integrand is given in Figure M. 



Using (|4.l| ) we have e.g., Koi{t,s) = K{t,ri{s),9i{s)](j)Q) where Vi and 9i are obtained 



using elementary trigonometry. We shall use the term "wedge kernel" to refer to both (4.1) 
and Kfj_y{t, s) as there is little danger of confusion. 




2-K + ia 



2-K — ia 



Figure 8: Wedge kernel: Analytic structure of integrand and integration contours. The inte- 
grand decays at infinity in the lightly shaded regions. Heavy shading denotes branch cuts in 
X{io) and crosses represent the poles of P{9,<j),uj). The dotted contour is the Carslaw contour 
C. For the semiclassical evaluation it is deformed to the solid contour which goes along the 
paths of steepest descent though zbvr. 
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4.1 How to calculate the contour integral 



As with all special functions, the calculation of (O) has to proceed through controlled ap- 
proximations. In general, a given approximation technique is only valid in a particular region 
of parameter space so that one must consider a range of techniques in order to be able to cal- 
culate the function globally. For (|4.1| ) we basically have two regimes: the semiclassical regime, 
\^/z{xf — Xj)| S> 1, in which the method of steepest descents is suitable; and the quantum 
regime, \^/z{xf — Xj)| <C 1, where the standard eigenfunction series is suitable. By developing 
each technique sufficiently one can extend each of them, with sufficient computational effort, 
so that the two regimes overlap. We are then able to calculate ( [4.1D for all parameters. 

Naturally, we are initially interested in an effective semiclassical approximation. Previous 
analyses [46, ^] were concerned with the scattering problem. The analysis for the Greens 



function is similar. However, because we require a controlled approximation it is neccessary to 
undertake a somewhat deeper treatment. Let us begin by outlining the standard semiclassical 
analysis. 

4.2 Asymptotic analysis and geometrical optics interpretation 

The semiclassical limit corresponds to — > co. In this limit the integrand is dominated by 
gi/^A(a;) gQ ^-^^^ ^j^g method of steepest descents is appropriate. One readily locates saddle 
points at zbvr together with the local path of steepest descent. Deforming the Carslaw contour 
C to the steepest descent contours one collects contributions from the real poles in the region 
— vr < Re(Lij) < vr (Figure ^). The pole contributions are of geometrical optics form, i.e. a free 
Green's function in which the distance, |xj — Xj|, is generalised to include paths from Xj to 
Xf via classically possible collisions with the walls. All possible classical paths are represented 
with an appropriate phase for the corresponding free Green's function. The remaining contour 
integrals, in the steepest descent approximation, can be interpreted in terms of a classical 
path from Xj to xj via the vertex. The importance of such paths, which are not classical 
but nevertheless minimise the action, was recognised by Keller ^ in his geometrical theory of 
diffraction. 

The steepest descent integrations are not uniformly of Gaussian form. The poles, located at 
2n(j) — 6i, (n e Z), pass through the the steepest descent contour at certain values of 6i so that 
the number of geometrical optics terms varies. This simply corresponds to the motion of the 
"shadow" of the source point. This apparently discontinuous behaviour is of course an artifact 
of the method of steepest descent. It is well known that the geometric terms "switch on" with 
an error function modulation. This is the Stokes' phenomenon and a uniform approach to it 
has been studied in recent years pO]. A uniform approach to a special case of the scattering 



problem was given by Pauli |46] 



It would be natural at this point to derive asymptotic series. This only requires a local 
knowledge of the integrand (in particular the steepest descent contour) at the saddle points 
and at infinity. In order to determine the nature of the error and, if possible, bound it, this 
understanding needs to be global. Because such an analysis is rarely presented we shall give it 
in some detail. 

While we could in principle carry out an error analysis for the asymptotic series, this is not 
really neccessary for our purposes. Asymptotic series are useful to understand the asymptotic 
form of the results. We already understand this through the standard interpretation in terms 
of (generalised) classical paths. All we require is well-controlled approximations in order to 
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validate the non-singular integral equation. This is conveniently and practically obtained via 
Gauss-Hermite quadrature along the steepest descent contour. Numerical quadrature along 
any other contour gives poor convergence due to oscillations in the integrand. (This is of course 
the rationale behind the steepest descent method.) We only know of one other study which 
uses the method of steepest descents as the basis of a numerical quadrature |48, p. 265]. 



4.3 Global understanding of integrand: u; ^ t transformation 

The steepest descent contour is defined implicitly by 

\{uj) = Ao H , -oo <t <oo (4.3) 

where Aq = Vcosha + 1. Rather than regard the variable t as a parameterisation of the 
steepest descent contour, consider ( [4. 3D to define a transformation to a complex i-plane. In 
the t-plane the contour is simply along the real axis. The integrand is of a semiclassical form 
only if it is a well peaked Gaussian in the t variable. 

Let us first write the conformal transformation (^]^) in the following dimensionless manner 



it^ L l + cos^_^ ^^^^^ 



a 

where 

r] = f^XQ = ^(n + rf) (4.5a) 

(7 = Aq = 1 + cosh a (4.5b) 

r] is essentially a dimensionless wave vector. The semiclassical limit is conventionally regarded 
as S> 1, corresponding to a de Broglie wavelength which is small compared to typical distances 
in the problem. The parameter o" is a purely geometrical, scale invariant parameter, a ^ 1 
indicates that that one of the distances rj,rj is much larger than the other, i.e. that one of 
them is "close" to the corner. As we shall see the semiclassical limit requires not only the 
conventional t] ^ 1 but also requires cr <C 1. 

w <-> t can be constructed via the sequence of transformations uj ^ u, u ^ v, and v ^ t 
where 

1 + cos LO / , X 

u = (4.6a) 

a 

— = (4.6b) 
V 

= v (4.6c) 

The choice of variable u is suggested by the fact that n = at the saddle point and turns out 
to be convenient. The other choices are more or less natural. 

u; <-> t is summarised in Figure |9[ The sequence of transformations leading to it is described 
in Figures 3(a), |l^, ^ and 12. The transformation is visualised by following features in the 



u plane through the transformation sequence. It is also convenient to follow the quadrants 
of the u plane. Arrows are used to denote how various branches should be glued together. 
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For example in Figure 9(a) we have a arrow labelled with a solid circle which crosses the line 
Im(u;) = vr, Re(Li;) > 0. In Figure |l^ we show the corresponding motion between branches of 
u{uj). 

Let us first consider lo <^ u. The images of the u quadrants in the lv plane are determined 
in the following manner. The images of the real and imaginary u axes are readily determined. 
They partition the w-plane. One can identify which region corresponds to which quadrant of 
the li-plane by using the asymptotic result u = 1 + cos to ~ 1 + IgT^R'^C'^) as Im(cj) ±00 
and continuity. The transformation is clearly multivalued. Single valued transformations 
(branches) uj{u) can be obtained by restricting the domain of uj. This can be done in an 
infinite number of ways so that we have the luxury of choosing the most convenient. For 
example we can partition the uj plane into the regions nvr < Re(w) < (n + l)7r, n G Z. Each 
region maps into the n-plane in a 1 1 manner. Alternatively could choose the partition 
2mT < Re{uj) < 2{n + 1)tt, lm{uj) > and 2mT < Re{uj) < 2{n + l)7r, Im(a;) < 0, n G Z. 
Any set of curves can be used in to partition the lo plane provided the regions obtained each 
cover the entire u plane once. The convenient partition for our purposes is the first choice. In 



Figures 9(a) and IC we show the transformation for the n = and n = 1 branches. It will 
become apparent in due course that we need to follow two branches at a time. For the moment 
the reader should note that the steepest descent contour through vr occupies both branches. 

n <-> V is readily constucted using a sequence of elementary transformations . The w-plane 
is given in Figure 11. An important feature to follow is the branch cut at n = -. In the 



I'-plane this becomes a vertical branch cut with the branch point at 

2i\fzri for rf > ri 



2iy/zrf for rj > rf i"^-"^) 

All that remains is v <^ t. Although this is just a square root, some thought is necessary. 
Up until now we have used the symbol with minimal explanation. In general this symbol 
denotes the standard square root, which has a branch cut along the negative real axis. Let us 
define another square root whose branch cut is along the negative imaginary axis. As this can 
be a point of confusion let us be very precise 

Let z = \z\e'\ with - - < < — . Then ^\ . = \z\^l^i^l'^, (4.8) 

We map the n = copy of the v-plane to the t-plane using v = V^|_joo and the n = 1 



copy of the u-plane to the t-plane using v = — y/z\_^^. We show the results in Figure 12 
where it is clear that we can paste the two regions together. The arrows indicate that this 
can be done in an analytic manner to give a 1 <-> 1 map of < Re(tj) < 2tt ^ (most of the) 
i-plane. We denote this as a;o(i) and call it the principal branch. 

Let us then summarise algebraically. We transform from oo to t using 

1 + cos UJ , . s 

u = (4.9a) 

a 

V = -ir/(-l + y/l-u) (4.9b) 

t={ 0<Re(a;)<7r ^^^^^ 

\- V^Uoo TT < Re(a;) < 27r 
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(a) tj-plane: Some features in lo <-* u are also shown (see Figure |lO[). The n = (n = 1) branch of this 
transformation corresponds to < Re(a;) < tv {n < Re(aj) < 2n). The bold dotted lines are images 
of the imaginary u axis. The bold solid lines are images of the real u axis. These lines partition the 
oi-plane into regions. The correspondance between these regions and the quadrants of the ii-plane is 
given by boxed numbers. 




(b) t- plane. Vc = i\fz(;ri + r/ — \r j — ri|) 



Figure 9: u; <-> t transformation: n = branch. Corresponding features in the two planes 
are labelled identically. Labelled arrows are used here and in other figures to indicate how 
branches are connected. 
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(a) n = branch (b) n = 1 branch 



Figure 10: oo ^ u ^ v ^ t transformation: u plane. Uq = ^ 




(a) n = branch (b) n = 1 branch 



Figure 11: uj <^ u <^ v <^ t transformation: v plane. Vc = i\[z{ri + rj — |rj — rj|). Analytic 
continuation into the darkly shaded regions corresponds to analytic continution through the 
cuts at 2n7r it ia, n G Z in the w-plane. 




(a) n = branch (b) n = 1 branch 



Figure 12: u u v t transformation: t plane 
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and from t to uj using 



r (4.10a) 



, 2 

IV 



u = l-(l + -] (4.10b) 



^ ^ . arccos(na - 1) -f < arg(t) < f ^ 

* 27r - arccos('ua - 1) < arg(t) < -f 



where we have used the principal branch of a; = arccos(2:). This maps the z~plane into the 
region < Re(u;) < vr. The other branches of the partition are obtained in a similar manner. 
In fact one can readily see that Wn{t) = 2mr -\- wo{t), n £ TL. The steepest descent contour 
occupying < Re(Lt;) < 27r (with orientation given in Figure ^ is mapped to the path in the 
t-plane along the real axis from oo to — oo. 

We can now see that it is not just the poles which spoil the Gaussian form of the inte- 
grand. When mindzrjl, \ zr\) <^ 1 (i.e. o" ^ 1) the cut at -^/^ collapses into the origin spoiling 
the Gaussian form of the integrand. This is the origin of the more subtle, non-semiclassical 
behaviour previously mentioned. At first sight it might appear that our formalism is always 
semiclassical, with its natural cutoff guaranteeing that + has a minimum value of macro- 
scopic magnitude. However because rj and are not so bounded from below there remain 
benign but nevertheless fully quantum contributions in the corners. 

4.4 Uniform calculation of integral 

Deforming C to the steepest descent contour and changing integration variables to t, we readily 
obtain 

Pole Terms- K^°\rf,ri,0i) + K^-^\rf,ri,ei) (4.11) 

where 

/oo 
P^^\t)fit)e-'"dt (4.12) 
-oo 

and 

pW(t) = P(0„0,c^„(t)) (4.13a) 
fit) = -^tHi'\r^ + it')e'' (4.13b) 

The integrations (4.12) are now in standard Gaussian form. 



For our purposes the most convenient way to proceed is to evaluate the integrals numerically 
using Gauss-Hermite quadrature. The quadrature converges well if f{t)P{t) is a smooth, 
slowly varying function. Because our construction guarantees that rj is reasonably large, this 
is true for f{t). P{t) however has poles and this condition doesn't apply when the poles are 



close to the real axis. This corresponds to Xf being in the penumbral region [49| so that the 



Gaussian approximation breaks down (see also [35, Figure 11.9]). We found that the presence 
of "close" poles spoiled the quadrature to an unacceptable degree for a significant number 
of matrix elements. It thus seems neccessary to handle such cases separately using some an 
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alternative numerical method. Fortunately, as we have mentioned, this is essentially a Stokes' 
phenomenon so that we are able to isolate the penumbral behaviour using well known uniform 



approaches pO, 51, p 



The basic idea is to analytically subtract the poles out of P^'^\t). This can be done a 
number of ways, the simplest being a recursive construction. Suppose we wish to subtract out 
poles at {llIq, ti^i, 1^25 • • • ,^n}- Let = t{ujn) and form (temporarily dropping the superscript) 

Q^{t) = P{t){t-to) (4.14) 

Defining 



Po = Qoito) (4.15a) 
'it) - . 
t-to 



Roit) = (4.15b) 



we obtain 



P{t) = + Roit) (4.16) 
t — to 

The pole at to is thus isolated as a simple pole in the first term. We can recursively subtract 
off other poles. Suppose that — 1 poles have already been subtracted off. Form 

QNit) = RN-imt-tN) (4.17) 

and define 

Pn = QNitN) (4.18a) 

Rj,{t) = (4.18b) 

we readily obtain 

N p 

Pit) = Y.-^ + RN{t) (4.19) 

n=0 " 

Although it is not immediately apparent from the construction, f„ is the residue of P(t) at 
tn- This is obvious from other arguments in complex analysis. Evaluating the residue we get 

„ u(j) sin uJn / , X 

It is natural to order the poles according to their proximity to the real axis and to subtract 
them out in this order. Plotting Rn{t) after each subtraction we noted the smoothing effect. 
As the proceedure is straightforward all the poles in the t-plane were subtracted out . 

We now proceed with the uniform analysis of the wedge kernel. Deform the Carslaw contour 
toward the steepest descent contour but this time don't cross the poles. This gives two contours 
To and T_i which are totally contained in the n = and n = — 1 branches of the w-plane. 
Rewrite (p^) as 



I^-^\rf,n,ei)-I^''\rf,n,ei) (4.21) 



28 



where 

I^^\rf,n,ei)= I P^^\t)f{t)e-''dt (4.22) 



The image of Tq (T_i) in the t-plane is a contour which goes from — cxd to oo passing over 
(under) all the poles under (over) the cut at (— )y^- 

Using ( [4.19| ) and performing one last one last pole subtraction we can write the integrand 
for n = as 

P^'Ht)f{t)= Yl E pW ^^^|:f(^"^ +<)(0 (4.23) 

0<Re(a;„)<27r 0<Re(aj„)<27r 

with a similar result for n = — 1. The integration is analogously decomposed into three terms. 
The absence of poles in the last two terms allows us to deform the contour Tj to the real axis. 
Using 



To t 



iiiw{z) z above T_i (4.24a) 
-iiTw{—z) z below Tq (4.24b) 



(Ref. @, Chapter 7) we finally obtain 

K{rf,rifii;4>) = 

0<Re(w,i)<27r -27r<Re(<^,i)<0 

0<lte^.)<2, (4.25) 

-27r<Rc(a;„)<0 "^-"^ " 

/CO . , 
(<-^^(t)-<°).,(t)) e-'^dt 
oo 

where A^^*^) (A^(~^)) is the number of poles in the n = (n = — 1) t-plane. 

(|4.25| ) is quite an interesting form. In the first two sums the contributions from the poles 
satisfying — vr < Re(LL'„) < tt are the geometrical optics terms, modulated by an error function. 
In this regime the asymptotic expansion of the error function has an exponential term which 
combines with the other factors to give the standard geometrical optics term. The first two 
sums also contain other poles. In the region containing these poles the exponential term is 
absent from the asymptotic expansion of the error function due to Stokes' phenomenon and 
the resulting terms in the sum are of higher order. We shall refer to the first two sums as the 
uniform geometric contribution to the wedge kernel. 

The "pole near a contour" problem is now isolated in the (well documented) error function. 
The remaining terms contain relatively smooth integrands and can be evaluated in a straight- 
forward manner using Gauss-Hermite quadrature. They can also be evaluated asymptotically 
by expanding the integrands around the origin. 
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A priori, one would now expect to have to consider the case min(|zrj| 



\zr\ 



<C 1. As is 



clear from Figure 9(b) , the semiclassical approach is not suitable in this "quantum" regime in 
which one should use the standard eigenfunction series. In Section ^ it will turn out that the 
parameter min(|2;rj|, \ zr\) never gets small enough to require the use of this series. 



5 Validation: Numerical Results 

5.1 Numerical discretisation 

The solution of integral equations via quadrature is a standard numerical technique described 
in many textbooks and with a vast research literature. The basic idea is to first choose some 
quadrature rule to approximate the integration. From here one readily obtains an inhomoge- 
nous linear equation for the unknown function at the quadrature nodes. This equation can 
then be solved using standard methods of numerical linear algebra. One then increases the 
number of nodes to achieve sufficient convergence. 

In this section our main purpose is to validate the wedge-based approach. To do this we 
content ourselves with the determination of the first dozen or so eigenvalues for some typical 
triangles. This can be done comfortably using a quadrature similar to that used by Berry and 



Wilkinson |27] and carrying out the calculation entirely within the Mathematica package [53]. 
Most cases could be done with an overnight run on a standard workstation. 

We expect that one should be able to use more sophisticated approaches very profitably in 
the wedge-based case. We don't expect this to be so for the singular case. 



5.2 Quadrature 

The geometric nature of the integration variable, being an arc length along the boundary, 
allows one to readily visualise the quadrature. In order to keep things simple, we consider 
cases where the boundary is a uniform 1-dimensional lattice. The length of each side is 
given as an integer multiple of some of a lattice spacing e and we give these lengths in the 
form of a triple, {no,ni,n2}. For example, the Pythagorean triangle with sides 1, | and 
I can be specified by the triple {3,4,5} with a lattice spacing of e = |. This defines 12 
quadrature nodes. We take {no,ni,n2} and e to define a base quadrature. We refine this 
quadrature to df{no,ni,n2} and e/df where the discretisation factor df is taken through the 
values 2,3,4, .. . in order to demonstrate convergence. For the wedge-based case we need to 
specify a subdivision point for each side. This is readily accomodated in an obvious manner 
with the notation {{no,ni}, {^-2,^3}, {714, 715}}. 

The quadrature nodes are chosen to lie at the midpoint of each unit cell so that the 
simplest quadrature is the midpoint rule. Our initial instinct was to avoid using the corners 
as quadrature points. 



5.3 Computational effort 

In order for the quadrature to be a good approximation it is necessary that one chooses enough 
nodal points so that the oscillations in the kernel are well sampled. The kernel oscillates with 
a wave number of the order ^/zL where L denotes a typical length scale for the triangle. A 
successful quadrature then requires the total number of nodes nj to be much larger than ^/zL. 
Li and Robnik ||3^ use a parameter b = 2Trnt/{y/zL) which gives the number of number of 
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nodes within one de Broglie wavelength. An effective discretisation requires rit ^ y/zL or 
6> 1. 

A discretisation of order nt requires 0{nf) operations to evaluate the kernel and then 0{n^) 
operations to evaluate the determinant. In our work the Bessel function evaluations were the 
most expensive part of the calculation so that effort grew quadratically in y/zL. As higher 
energies are considered the effort will eventually cross over to the cubic regime. 



5.4 Results 

The eigenvalues are determined from the zeros of the Fredholm determinant, det(l — i^), which 
we approximate using the discrete form of the kernel defined by our quadrature. 

We consider in detail two cases: the equilateral triangle and the Pythagorean 3-4-5 triangle. 
The former case is an integrable billiard and the eigenvalues are known exactly. The integrable 
cases are described in [^]. The latter case was chosen as a "typical" non~integrable billiard. 

Let us first give an overview of the lower part of the spectrum for each case, presenting 
plots of the modulus of the Fredholm determinant for real energies. The cases considered are 
given m Table |l[ 



Abbreviation 


Base Quadrature 


free-equilateral 


{4,4,4},e = i 


free- 345 


{3,4,5},e = i 


wedge-equilateral 


{{2,2},{2,2},{2,2}},e = i 


wedge- 345 


{{2,l},{2,2},{3,2}},e = i 



Table 1: Base quadratures for cases considered together with a shorthand description. 



In all plots we have applied an appropriate energy scaling. For the integrable cases we 
have used a scaling which allows easy comparison with exact results. For the non-integrable 
cases we have used Weyl's result for the density of states to scale the energy so as to give, on 
average, a unit eigenvalue density. 

In Figure 13(a) we give the results for the free-equilateral case. The first thing one notices 
is that the determinant does not converge with increasing b i.e. the continuum limit does not 
exist. Despite this, the variation of the determinant still clearly resolves the eigenvalues. The 
spectrum, given in appropriate units hy E = m? + — mn, 1 < n < m, is reproduced. 
Cusp behaviour near the zeros indicates isolated eigenvalues and smooth behaviour indicates 
degeneracies. 

The wedge-equilateral case is shown in Figure 13(b). Here one observes the convergence 
that Fredholm theory demands. One expects that the lack of convergence in the free-based 
case is due to the singular nature of the integral equations. 

In Figures 14(a)| and 14(b) we show the corresponding results for the 3-4-5 triangle. Once 
again the determinant doesn't converge in the free case. We have used the uniform-geometric 
component of the wedge kernel. As we shall see this gives accurate results. 

In order to examine convergence more closely in the free-based approach, we show in 
Figure 15 a log — log plot of the absolute value of the determinant versus the discretisation 
factor at some particular values of the energy for the 3-4-5 triangle. It is clear that the 
determinant decreases with df in a power law manner with a exponent ~ —1.5. The exponent 
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(a) Free-based Predholm determinant 




3 7 1213 19 21 

Scaled Energy, Eg 

(b) Wedge-based Fredholm determinant 



Figure 13: Equilateral Triangle: solution of singular and non-singular integral equations by 
numerical discretisation. The energy has been scaled to facilitate easy comparision with the 
exact solution. The number of lattice spacings per deBroglie wavelength, b = 6df/^/E^. 
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(a) Free-based Predholm determinant 
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(b) Wedge-based EYedholm determinant: uniform geometric kernel 



Figure 14: 3-4-5 Triangle: solution of singular and non-singular integral equations by numerical 
discretisation. The energy has been scaled to given a unit density of states. The number of 
lattice spacings per deBroglie wavelength, b = 4.3df/y 
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Figure 15: Convergence in free-based case: log — log plot of the modulus of the Fredholm 
determinant as a function of discretisation factor df, \ det{df)\, at fixed energy for the 3-4-5 
triangle. 



appears to be relatively uniform in energy. One could thus arrest the lack of convergence with 
an ad-hoc "renormalisation" factor. 

In Figure ^ we quantify the convergence of the wedge-based approach by showing a 
log — log plot of the first difference of the discretised Fredholm determinant versus dj for 
the equilateral triangle. The power law exponent is approximatly -3, which is consistent with 
the use of the midpoint rule. The straightness of the plot indicates that standard methods to 
accelerate convergence should work well. 

In Figure 17 we give detailed results in the vicinity of the 3rd eigenvalue of 3-4-5 triangle. 
We compare the geometric, uniform-geometric and full kernel-based approximations for the 
wedge-based method together with 3 discretisations of the free-based result. We see that the 
approximate wedge kernels give acceptable accuracy. We also note that the wedge-based and 
free-based methods agree very well in the location of the zeros. 

To obtain 6 figure accuracy uniformly across all matrix elements for the full kernel we 
used an 8-point quadrature for the Gaussian integrals in ( 4.25| ). This is sufficient given the 
discretisation error in the numerical solution of the integral equation. This quadrature is 
relatively expensive due to the need for extra Bessel function evaluations. 

As we mentioned at the end of Section the semiclassical form of the wedge kernel is all we 
need. Because of the midpoint quadrature there is a minimum value of rj and r^. Our modest 
numerical goals meant that the discretisation remained fairly coarse so that min(|zrj|, |zr|) 
was never small enough to see the quantum regime. 

In all figures one observes an expected decrease in the quality of eigenvalue determination 
as the quality of the discretisation b decreases. Where direct comparision is possible the 
wedge-based code seems to resolve the eigenvalues better for a given b. 

In some of the figures the curves are not perfectly smooth. This is due to the limited 
resolution of the plot. Each plot displays the determinant sampled at about 100 equally 
spaced points. The samples are connected with lines by the plotting routine gnuplot. Where 
the determinant varies quickly, such as near an eigenvalue, a higher resolution is needed. It 
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Figure 16: Convergence in wedge-based case: log — log plot of the first difference of the Fred- 
holm determinant as a function of discretisation factor df, \det {df + 1) — det {df)\, at fixed 
energy for the 3-4-5 triangle. 



is best to offer a separate plot in such regimes, as we have done in Figure 0. Such plots are 



needed to see the actual degree to which eigenvalues are resolved. For example in Figure 14(b) 



the third eigenvalue does not appear to be resolved very well. This is misleading because the 
the only energies sampled in this vicinity are Eg = 5.267 and 5.333. Figure 17 shows the true 
story. Also the power law decrease of the free determinant with df artificially smoothens plots. 
This might falsely lead one to think that the free-based code gives better resolution that the 
wedge-based code. A "renormalisation" of the free results would give a fairer comparision. 



5.5 Numerical checks 

In order to develop confidence in our results we applied several checks in order to diminish the 
possiblity that numerical or coding errors were present. The first and most obvious check is 
the reproduction the exact spectrum for the equilateral triangle in both versions of the code. 
For the 3-4-5 case, and for nonintegrable triangles in general, there are no exact results for 
the eigenvalues. In these cases we apply consistency checks. Our main consistency check was 
to compare the results of the wedge-based code with the free-based code. We found that for 
the common eigenvalues, the agreement was very good. This can be seen in the case of the 
3th eigenvalue in Figure 0. As a further check we also considered the free-354 and free- 
534 triangles, obtaining results identical to the free-345 case. To eliminate the presence of 
integrable wedges, we also considered the 2-3-4 triangle. 

The lack of determinant convergence for the free-based code was originally thought to be 
possibly a result of coding error. The free-based code was originally intended as a "warm-up" 
exercise. After it was written we proceeded onto the more involved wedge-based calculation 
and improved the general quality of our software considerably. It was then relatively easy to 
modify the wedge-based code to produce, independently, a new version of the free-based code, 
which reproduced the original results. We further tested the free-based code by considering 
the other two integrable triangles: the isoceles right angled triangle (45,45,90) and the half- 
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Figure 17: Fredholm determinants in the vicinity of the 3rd eigenvalue of the 3-4-5 triangle: 
detailed comparision of singular and non-singular approaches. For the latter the effect of 
approximations in the wedge kernel is also shown. F = free-based; WG = wedge-based, 
geometric kernel; WU = wedge-based, uniform geometric kernel; W = wedge-based with full 
kernel. 



equilateral triangle (30,60,90). Our present code is written for triangles with rational sides. 
We thus considered the {12, 12, 17} and {15, 26, 30} triangles, with suitable values of e. These 
approximate the two cases to 1 part in 100. We found full agreement with the exact spectra. 



6 Periodic Orbit Expansions 

We now have a mathematically well behaved reformulation of the boundary value problem 
which has been well validated using numerical techniques. 

An intuitive understanding of the high energy, semiclassical eigenvalues could now be ob- 
tained by deriving a periodic orbit expansion using the Fredholm series, much as Georgeot and 
Prange describe |2^ ]. The calculation of the periodic orbit expansions is now a matter of using 
the semiclassical form of the wedge kernel to calculate Tr„ asymptotically. Similar analysis has 



been done in other contexts |29, ^]. In our case the analysis will differ due to the subdivision 
of the sides of the polygon and the uniform approach to the diffractive contribution. 

It is more important to first examine the nature of the periodic orbit expansions using 
Fredholm theory. In particular we wish to compare the computational effort of the periodic 
orbit expansions with standard numerical discretisation. 

6.1 Convergence of Fredholm series 

Standard Fredholm theory |23, §6.5] shows that series for the Fredholm determinant, 

oo 

A{z) =J2 ^n{z), (6.1) 

n=0 
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is convergent by bounding |(5„| using Hadamard's inequality. From this bound one can readily 
estimate that \6n\ increases until 

Uc^WKf (6.2) 

after which it decreases. In practice one might then expect that at least Sric terms are needed 
in the Fredholm expansion. 

In order to apply this result we need to understand the behaviour of the norm for all values 
of the energy. In particular it is clear that one has an effective semiclassical approximation if 
and only if the norm decreases with energy. The faster the rate of decrease the better. 

It is relatively easy to estimate the norm in the semiclassical limit ^/z\A.^^\ ^ 1. In this 
limit the wedge kernel is dominanted by geometric optics terms so that it has the form 

i.„„(M)~-^.". 2 ^ (8.3) 

classical 1^^' ^ 
paths a 

The number of classical paths varies with t and s as the shadow of the source point moves. 
Diffractive contributions decay with a faster power law 
The semiclassical norm then becomes 

^y.y. 1 ds 1^'^ dt\K^^{t, s)|2 (6.4a) 



(6.4b) 




A" 11/21^/3 11/2 |A" I I 

= \fzL + oscillatory terms (6.4c) 

where L is length of order the size of the triangle. The oscillatory terms come from a steepest 
descents evaluation of the integration over the second sum. 

(]6.3|) is not uniformly valid for all t and s. Although our technique avoids the unit source 
singularity there remains non-classical behaviour when one of \fzri or \fzrf is very small. 
One expects that (|6.3| ) will cross over into some other form in this boundary layer region. We 
have no reason to believe that incorporating the boundary layer contribution will dramatically 
alter the above estimation of the norm. To confirm this we generated log — log plots of the 
norm versus energy and confirmed the above behaviour. However our discretisation was not 
fine enough, nor was our energy high enough, to see any boundary layer behaviour in the 
calculation of the kernel. It is thus conceivable that at higher energies some change in the 
norm behaviour may be observed. 

We thus have 

nc ~ \fzL (6.5) 

so that the number of terms needed grows with the square root of the energy. This accords with 
the literature on periodic orbit expansions but is not what one would expect of a semiclassical 
expansion. 
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6.2 The correspondance principle 

In the old quantum theory Bohr presumed that any semiclassical approximation should improve 
in "the limit of large quantum numbers" . In our context one would expect that the semiclassical 
approximation should improve as energy increased. This is clearly not the case for periodic 
orbit expansions. The semiclassical Fredholm series is useful only if ^/zL is large enough to 
use the semiclassical approximation for bn{z) but small enough so that not too many terms 
are needed. 

We can offer an elementary analogy to impress apon the reader the true nature of periodic 
orbit expansions. Consider the asymptotic calculation of the zeros of the Bessel function, i.e. 
solve Jv{,z) = for z ^ oo. The Bessel function can be defined in terms of series or integral 
representations. It is natural to handle this problem using the standard asymptotic result 

Pi 1 1 

Jy{z) ^ \\ — CO^{z- -VTX - -TX) (6.6) 

V vrz 2 4 

Zeros obtained in this way are very accurate and improve in accuracy as z — > cxd. In actual 
fact the accuracy is acceptable right down to the ground state! ([^, 9.5.12]) A less natural 
way to proceed would be to use the series representation 

oo / \ 2\k 

because the number of terms need to obtain good convergence is of the order of |z. Finding 
eigenvalues from periodic orbit expansions is like using the power series expansions to find 
zeros of Bessel functions. It is clear then that 

The problem of semiclassical quantisation of chaotic systems will not be solved until 
one has found the analog of (|6.6| ) . 

6.3 Computational effort: comparison with discretisation 

Let us now examine the computational effort needed to evaluate the semiclassical Fredholm 
series. 

The number of terms need to obtain good convergence is of order the dimensionless 
wavenumber \/zL. A semiclassical calculation of 5„ involves all periodic orbits with n bounces 
for the triangular billiard. The number of such orbits is known to grow polynomially in n [2 5( 1. 
The order of this polynomial is not known but there are some conjectures. For integrable tri- 
angles it is relatively easy to see that the growth is linear. Either way the computational effort 
required to calculate high energy eigenvalues grows polynomially in ^/zL. Discretisation of the 
boundary integral equation requires computational effort which grows like {^/zL)^ although in 
our case, as we discussed in Section ^, it was effectively quadratic. This becomes comparable to 
the semiclassical Fredholm series for integrable billiards. For generic billiards we would expect 
the effort to be qualitatively comparable 

This is in contrast to the case of chaotic smooth billiards. Our estimate of the norm in the 
semiclassical limit should essentially carry over to this case. Following Georgeot and Prange 
it is clear that exponential proliferation causes the calculation of semiclassical Fredholm series 
to grow exponentially. It is then clear that for smooth billiards the boundary integral method 
is much more efficient. It thus appears that 
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Even with convergence questions settled, the problem of exponential proliferation 
remains. 

We shall discuss this point further in Section |^. 

7 Discussion 

7.1 Singular integral equation 

7.1.1 Analogy with quantum field theory 

The divergence of the norm in the free-based integral equation suggests that the Fredholm 
series itself may contain divergences. In fact, using an argument virtually identical to the one 
we used for the norm, one can show that Tr{K'^) diverges logarithmically. One can extend 
this argument to show that the traces of all even powers of K contain logarithmic divergences. 
These divergences are contained in contributions which scatter back and forth between two 
given sides. It appears there are no other more serious divergences. There is an obvious 
similarity between this behaviour and that of perturbation series in quantum field theory. 

As in field theory, these divergences are due to the short distance behaviour of the theory. 
The semiclassical contribution comes from the long distance behaviour y/zL ^ 1. Because 
these two regimes are well separated it is easy to overlook the divergences in a semiclassical 
calculation. However, we know from the literature, in particular the work of Pavloff and 
Schmitt [^], that periodic orbit expansions do give a reasonable degree of agreement. This 
suggests that, like in field theory, the semiclassical properties may be insensitive to the short 
distance divergences. Because of this we were tempted, at various stages of our study, to 
consider some sort of "renormalised" periodic orbit expansion for the free-based approach. The 
removal of singular behaviour in perturbation series is usually done in an ad-hoc manner by 
"regularising" the divergences with artificially added convergence factors. These depend on a 
cutoff parameter 1 /A in such a way that A — > oo recovers the original integrand. One evaluates 
the integrals for finite A, artificially removes the singular terms by some "renormalisation" 
proceedure, and then takes A — > oo to obtain "finite" results. 

The natural regularisation for our problem is to introduce a short distance cutoff, much as 
we do for the estimation of the norm with the definition of the annular region. Alternatively we 
could apply a "dimensional regularisation" by replacing the H^^^ in the free-Green's function 
with hP. 

One then needs to develop a rationale to treat the divergent parts of the traces, that is a 
"renormalisation" process. There is nowhere obvious to "absorb" the divergences. Furthermore 
it is not clear how this process would modify the periodic orbit expansion. Rather than proceed 
along this line we developed the wedge-based approach. 

We have no doubt that the pathologies in the singular formulation have physical interpre- 
tations. In particular we suspect that other boundary value problems are mixed in as a result 
of the the sloppy treatment of the corners. Perhaps this leads to a leaking of "flux". 

7.1.2 Wedge based approach as renormalisation 

Our wedge-based construction clearly has the spirit of a "regularisation" because the division 
of each side into segments is clearly some sort of "cutoff". However our approach is not a 
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post-facto arbitrary alteration of the result, it is a first principles ab-initio construction. 

The need for a wedge-based construction is also indicated by the nature of the divergences 
in the singular approach. As we indicated the divergences came from repeated propagation 
between two given sides. In the perturbation series culture one would have removed these 
divergences by summing over "ladder" diagrams. We can interpret the perturbation about 
the wedge solution as a way of summing over these diagrams. The Faddeev equations have a 
similar interpretation. 

Even though the wedge-based approach is clearly preferable there are circumstances where 
such an approach is not available. Our study provides a simple model in which one can critically 
examine the concept of "renormalisation" in more depth. 

7.1.3 Wedge based approach and renormalisation group 

When we first came across the need to subdivide the sides of the triangle, we were quite 
uncomfortable about the arbitrariness of the subdivision. We have come to think that this 
subdivision cannot be avoided. In actual fact there is a natural choice of subdivision: given the 
role of the norm in controlling the convergence of the Fredholm series it is natural to minimise 
the norm with respect to this subdivision. 

In many ways this is like the field theoretic renormalisation group program. Here one 
introduces an arbitrary parameter, the renormalisation point /x. This is where one fixes coupling 
constants. Physical quantities should be independent of /x and this is imposed in the theory 
by setting the first derivative with respect to fi to zero. In our case setting the derivative to 
zero is clearly associated with improving convergence. 

7.2 Smooth billiards 

7.2.1 Doubts concerning the standard approach 

The work of Georgeot and Prange suggests that one can immediately obtain semiclassical 
Fredholm series for smooth billiards from the integral equations of Balian and Bloch. Our 
experience with the triangular billiard has led us to entertain doubt and conjecture that these 
non-singular integral equations are incorrect. As we hope to have impressed upon the reader, 
the major conceptual difficulty in our work has been the determination of mathematically well 
behaved integral equations. 

We have two pieces of evidence for this conjecture. First of all, a careful boundary limit 
along the lines of Section § seems to produce e dependence similar to that in the corner 
boundary layers in the triangle. Secondly, there are problems in the standard approach for 
non-convex domains. Li and Robnik ||3^ have shown that the boundary integral method fails. 
Balian and Bloch have shown that classically forbidden chordal paths in the periodic orbit 
expansion cancel II]. In contrast there is no difficultly with non-convex polygonal domains 
in our approach. Some matrix elements will have no geometric contributions, but there will 
still be diffractive contributions which creep around the corners much as one would expect. 
We thus feel that the kernel for smooth billiards should have a diffractive component in order 
to accomodate the non-convex case. 

We suspect the Balian and Bloch calculation is correct only for sufficiently smoothed density 
of states expansions. In fact this is all they claim, but their work seems to be regarded in the 
same light as the Gutzwiller trace formula, i.e. capable of determining individual eigenvalues. 



40 



We thus feel that a critical look at smooth billiards is needed. In particular we offer some 
alternative approaches. 



7.2.2 Alternative approaches to smooth billiards 

With suitable solutions to perturb about, our construction can be readily generalised to chaotic 
billiards. Given a suitable form of the Green's function outside a disk, one could readily obtain 
non-singular integral equations for the 3~disk system. Here the finite separation of the disks 
gives a lower bound on propagation distances so that the resulting integral equations will 
clearly be non-singular. Given the Green's function for the hyperbola wedge, which might be 
derived from eigenfunction expansions, we could obtain convergent periodic orbit expansions 
for the hyperbola billiard. 

An interesting and more general proceedure would be to approach smooth billiards using a 
suitable sequence of polygonal approximations. One could take such a limit before or after the 
integral equations are solved. A natural starting point would be to consider the regular n-gon 
and take the limit n — > oo to see what sort of integral equation one obtains for the circle. This 
approach has been used by Kac [^] to derive the curvature dependent constant term in the 
Weyl formula for smooth billiards. Before using this approch he tried to "perturb" about the 
solution of the diffusion equation on a disk. He was forced into using a polygonal approach 
due to the lack of a "workable" solution to this problem. 



7.3 Modified periodic orbit expansions 

Whether or not the standard integral equations for smooth billiards need to be replaced, much 
of our discussion is relevant to more general systems. Let us compare then our approach with 
the various modifications of the original Gutzwillers trace formula. 



7.3.1 Surface of section techniques 

The surface of section approach of Bogomolny |0, 21 1 takes a "section" in configuration 
space and expresses the eigenvalues as the zeros of a Fredholm determinant det(/ — T). T is 
an operator whose co-ordinate space representation is the semiclassical propagator containing 
all classical paths starting and ending, but not otherwise crossing, the surface of section. This 
method is only partially semiclassical since the Fredholm determinant has to be determined 
by a quantum technique: either discretision or diagonalisation in some basis. Bogomolny's 
technique is generally applicable and is not restricted to billiard problems. It was, however, 
motivated by the boundary integral approach to billiard problems. 

Our approach is a semi-rigorous realisation of Bogomolny's ideas. The semiclassical ex- 
pansion of the wedge propagator is clearly analogous to the T matrix. Bogolmony's approach 
corresponds to treating the wedge kernel semiclassically while treating the Fredholm deter- 
minant fully quantum mechanically. In our approach the Fredholm determinant arises as a 
natural consequence of the reformulation of the problem as an integral equation. 



As we discussed in Section 3.3 the surface of section/boundary integral equation method 
requires much less computation effort. The sacrifice of classical intuition in treating the Fred- 
holm determinant quantum mechanically seems a necessary price to pay in order to avoid 
exponential proliferation. 



41 



While most of the hterature views the problem of quantum chaos as determining quantum 
mechanical quantities from classical ones, some authors also considers the converse: studying 
classical mechanics using quantum mechanics. This point of view arises from regarding the 
semiclassical Fredholm series as a generating function for classical periodic orbits. One might 
regard quantum mechanics as a rival of the thermodynamic formalism, which considers similar 
generating functions. One might even muse that if experiment hadn't demanded it, one might 
have had to develop quantum mechanics as a linearisation of classical mechanics. From this 
point of view it appears that Bogomolny's hybrid method may be the most effective way to 
study both quantum and classical mechanics 



7.3.2 Cycle expansions 

As Georgeot and Prange have pointed out, the pseudo-orbit expansion, which arises from the 
expansion of the zeta function, is the Fredholm series. In a similar manner cycle expansions 
are also the result of the expansion of zeta functions. The way in which terms are grouped 
together so that longer orbits are grouped with shorter ones is very similar to the form of the 
Fredholm series I^Sl. 



Cycle expansions, however, have a further rationale within classical mechanics [56|. The 
basic element of cycle expansions is that the dynamics is topologically organised in terms of 
a symbolic dynamics. This involves a considered partition of the phase space. Each region is 
labelled with a symbol so that each orbit, in particular each periodic orbit, has a symbol se- 
quence associated with it. This topological understanding allows one to construct well behaved 
cycle expansions of the Fredholm determinant. The cycle expansions themselves come from 
expansions of determinants associated with automata describing the grammar of the allowed 
symbols. Our wedge construction is similar in spirit in that the phase space (the boundary) is 
divided into regions. In our case the form of the division is not dictated by the classical dynam- 
ics but by the need to generate non-singular integral equations. We feel that there should be a 
relationship between these different criteria. In particular we feel that the quantum mechanics 
of a problem might provide a useful guide to developing suitable symbolic dynamics for the 
classical problem. 



7.4 Precursors of the Fredholm series 

Whilst the use of the Fredholm theory began with Georgeot and Prange, the Fredholm series, 
as they point out, has already appeared in the guise of the pseudo-orbit and cycle expansions. 
It has thus been a part of the quantum chaos literature for some time. As another example 



Alsonso and Gaspard use the logdet = Trlog identity |38, Eqn. (2.8)] in a Balian-Bloch 



context. Bogolmony |55, Eqns (23) and (41)] gives arguments for the convergence of the 
Fredholm series in the context of his surface of section technique without specifically refering 
to Fredholm theory. 

In much of the literature emphasis is given to the zeta function, probably as a result of 
the infiuence of the number theoretical and hyperbolic geometry literature. This influence 
then motivates the "solution" of the convergence problem by formally expanding out the zeta 
function and grouping the terms in a well motivated, but essentially ad-hoc way. The Fredholm 
theory, with its emphasis on the integral equation reformulation of the boundary value problem, 
provides a direct and mathematically firm way to circumvent this roundabout procedure. 
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7.5 Improving convergence? 

Exponential proliferation is only a problem in the semiclassical Fredholm series because the 
norm of the kernel grows like \fzL. The question arises as to whether it is possible to further 
manipulate the integral equations to reduce this norm. In the literature on integral equa- 
tions it sometimes happens that iterating the integral equations can help. This leads to the 
consideration of the Fredholm Determinant det(l + IT"). It is fairly clear, however, that the 
semiclassical form of -fC" is of similar character to that of K so that the argument we presented 
for the norm should still apply. While we are not hopeful that such manipulations will lead to 
a norm which decreases with energy we feel that it is best to defer a full semiclassical analysis 
until this issue is fully investigated. 

7.6 Deriving Keller's geometrical theory of diffraction 

Our wedge-based approach is an example of a way in which one may place Keller's geometrical 
theory of diffraction Q on a mathematically firmer basis. Keller's theory begins with the stan- 
dard semiclassical solutions of the wave equation Ae^^ in which the action S is determined via 
an extended Fermat principle, which includes diffracted paths. Because these are approximate 
solutions to a dijferential equation they contain arbitrary constants. These have to be deter- 
mined using the boundary conditions. Keller does this by matching the approximate solutions 
to exactly known solutions near points of diffraction. The idea here is that in the short wave- 
length limit only the local geometry can effect the solution. Often the local geometry is one in 
which the problem has been exactly solved. One then expects to be able to obtain the local 
behaviour of the solution as a perturbation of this exact result. Matching both solutions to 
leading order will give the unknown constants. The spirit of Keller's technique is thus similar to 
that presented here. However, because we use Green's theorem and convert the problem to an 
integral equation, boundary conditions are automatically incorporated so that the somewhat 
awkward matching proceedure is not needed. Furthermore ours is an exact reformulation of 
the problem and the semiclassical approximation is done after the exact solution is obtained. 
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A Wedge Greens Function 

A statement of the boundary value problem for the wedge Green's function is obtained from 



(LI) by simply replacing the domain A with the wedge W = {{r,9) :0<6<((>,0<r< 00} 
and requiring that the solution decay at infinity. 

Carslaw |44] derived this Green's function using a method due to Sommerfeld [^]. The 



solution is constructed by altering a trivial contour integral representation for a suitable elemen- 
tary solution. While the motivation for this proceedure is somewhat obscure it is reasonably 
easy to verify the result. Here we offer an outline of the derivation which combines Carslaw's 



discussion with a method due to Williams |58]. Once again we shall try to maintain reasonable 
standards of rigor, but our discussion will not be rigorous. The reader can find yet another 
derivation in |47|. 
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We begin following Williams with the free Greens function in polar co-ordinates, 

Go{r,rf,e-ei) = -'-H^'^ (^VT ^r^ + rj-2rncos{e-ei)^ (A.l) 
We rewrite the Helmholtz equation which it satisfies in the following form 

ld_ / aGo\ ^ ^ (A 2) 

r dr \ dr J duj"^ 

where we have set 9 — Oi = lo, promoting it to a complex variable, 

One then writes down a angular superposition of elementary solutions as an ansatz 

Gw{r,e)= [ f{uj,e)Go{r,n-u:)duj (A.3) 
Jc 

The contour C is a suitably chosen path in the complex uj plane. Note the similarity to the 
method of separation of variables. 

In order to make practical use of this ansatz we must be able to interchange integration, 
differentiation and limits. Such manipulations are valid only if integral is uniformly convergent. 
The contour C must then avoid singularities of the integrand by a distance which is bounded 
below. Furthermore in any passage of the contour to infinity the integrand must go to zero 
sufficiently fast. 

We must also specify the branches of the square roots appearing in ( [A.l| ). This requires 
some anticipation of the ultimate form of the solution. From integrable cases and semiclassical 



arguments we expect to find contributions of the form ( |A.1| ) where Oi is replaced by angles 
reflecting image contributions. Here the second square root represents a distance and is thus 
real and positive. In order to satisfy the boundary conditions at infinity one then chooses a 
branch of the first square root to satisfy Im ^/z > 0. This has a branch cut along the positive 
real axis. 

The "distance" square root can be written as ^J2rir\{Lv). In the notation of section ^, 
the argument of the square root is essentially 1 — u. Using our knowledge of u; <-> u we can 
determine the analytic properties in the uj plane corresponding to different choices of branch. 

As our studies are limited to real energies our first choice was the branch Im-^/i > so 
that Go(r, rj;a;) decays exponentially in all directions as |u;| oo. Although this gives more 
choice for the contour C, it creates a branch cut in the w-plane starting at ia, going towards 
the origin along the imaginary axis, along the positive real axis until vr and then up the line 
TT + it, < t < OO. The other branch cuts are obtained by suitable refiection and translation. 

The image contributions, however, come from poles on the real axis so that this choice 
severly restricts (or at least complicates) the various contour deformations we need. To avoid 
this we choosing the standard square root, which has a branch cut along the negative real 
axis, to obtain a branch cut in the to plane which starts at ia and goes to ioo, again with 
others related by symmetry operations. The are depicted in Figure ^ together with shaded 
regions depicting where the integrand decays at infinity. We determine the latter by noting 
that ImA(u;) > requires 1 — to be in the upper half plane. Prom Figure |lO| one can see 



that this corresponds to regions |_3j and |_4j which we can trace back to the uj plane using 
Figure 9(a) . 
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Substituting the ansatz into the polar form of the Helmholtz equation, exchanging differen- 
tiation and integration, using Eqn. [A.2| to ehminate the r derivatives and integrating by parts 
twice, one obtains 



JC 



Go{r, ri;uj) doj 
B 



{r,ri;uj) - Go{r,ri;uj) — {u;,( 



(A.4) 



If C is chosen so that the latter "complex flux" term vanishes we see that the problem is 
reduced to a complexified 1-dimensional wave equation for f{ui,6) whose general solution is 
elementary 



/(tj, 9) = M{uj + 6) + N{uj - 6) 

Because the variables r and 9 are separated, the Dirichlet boundary conditions for 9 
can be satisfied with the choice 



f{uj,9) = M{uj + 9)-M{uj 
M{uj) = M(a; + 2(/)) 



(A.5) 
= 0,(/. 

(A.6a) 
(A.6b) 



The remaining function M{z) and the contour C must now be determined from the unit source 
condition. 

The unit source condition is satisfied if 



Gwir,9) = Go{r,9) + Gi{r,9) 



(A.7) 



where Gi{r,9) has no sources in the domain W. This can be arranged with a construction 
which is essentially the method of images. 
If 

1 



1 



2TTi Z + I 



(A. 



A.7 is satisfied if C can be deformed to 



hm M{z) = 

then M{uj — 9) has a pole of residue \ a,t uj = 9 
pick up this pole while ensuring that the remaining contour integral remains sourceless in W . 

In contemplating an analytic function M{z) satisfying Eqns. A.6b| and A. 8 we are naturally 
lead to consider the function 



M{z) 



lim 

Af^oo 2111 



1 ^ 

- E 



1 



2 + 



(A.9) 



It is possible to choose other functions by adding further poles in a periodic way. Such extra 
poles will generate features in the final solution which either spoil the uniform convergence 
or don't satisfy the boundary conditions. The principal value-type construction is needed to 
obtain convergence. The reader should note the resulting contour integral as a simple example 
of a conditionally convergent semiclassical series. This sum is a textbook example ([^, p.261]). 

^ 1 

lim > = 7rcot7r2; (A. 10) 



-N 



n 
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M{uj — 6) generates images at 9i — 2n(f). M{uj + 0) generates images at —Oi + 2n(/>. Some 
of these poles, with suitable contours, could generate unwanted source points images in the 
domain W. This can be prevented in a rather elegant manner by an appropriate choice of the 
contour C. 

We now have sufficient information to uniquely determine C. Requirements of uniform 
convergence dictate that C must not be near singularities for any value of 9i . In particular the 
contour cannot cut the real axis. The "complex flux" can then only vanish if the contour starts 
and ends at oo in the lightly shaded regions in Figure ^. We need loops in both the upper 
and lower half plane so that the contour can be deformed to the "direct" pole at to = 6 — 6i 
so that the unit source condition is satisfied. Furthermore the deformation must not pick up 
the poles in M{z) which correspond to unwanted source terms. A loop starting at Imw = oo 
in — 27r < Re a; < — vr and finishing at Ima; = oo in Stt < Re a; < 47r, together with it's lower 
half plane analog, for instance, will pick up unphysical sources. The only choice satisfying all 
these conditions is the Carslaw contour. 

The wedge Green's Function is thus 

Gw{r,9) = ^J^[P{6i-9;io)-P{9i + 9-io)]Go{r,rf,uj)dco (A.ll) 

where 



jTr{a+uj)/4> 



1 



gi7r(a+a))/</i ^ 1 f,—iTT{a+uj)/ip 



(A.12) 



The solution must satisfy reciprocity i.e. the solution should be symmetric on interchanging 
(r, 9) and (r^, 9i). This is not manifest but can be demonstrated by (i) isolating the part of the 
integral involving P{9i — 9;uj); (ii) deforming its contour C to one which is symmetric under 
Lo — > —uj; (iii) letting m — > —uj (C ^ C,Go Go^duj — > —du>) and (iv) using the result 
-P{a; -uj) = -1 + P{-a; uj). 

The wedge solution also gives us a deeper understanding of why the method of images 
is exact iff = ^, n = 1, 2, . . . . For these angles the integrand is 2tt periodic. The contour 
can then be deformed into two pieces which have opposite sense, related by a 2tt translation, 
picking up image points along the way. Periodicity causes the integrals to cancel. The classical 
problem is also integrable at these angles. 

A.l Wedge Kernel 

The integral kernels we want can now be calculated by differentiation. Any kernel ending on 
the lower edge can be expressed in terms of 



K{rf,0;ri,9,) = -VGw{rf,0;ri,9i) 

_ l_ dGw 

rj 89 



-{rf,9;ri,9i 



(A.13) 



Using 



lpi9,.) 



(A.14) 
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and integrating by parts, we get 



1 



<t>rf Jc 

The kernel ending on the upper edge is simply 
K{rf,(l);ri,9i) 



-VGwirf,(l);ri,6i) 

1_ dGw , 
V/ do 



(A.15) 

(A.16) 
(A.17) 



We proceed as for the other kernel, except that we need to use periodicity and change variables 
LU — >■ —CO in a manner similar to the reciprocity proof, to obtain 



4>rf Jc 4 



(A.18) 



This is identical ( A.15| ) except that the angle is measured from the ray = (f) m a, clockwise 
sense. This is in fact required by reflection symmetry 
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